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Mean spherical approximation (MSA), second-order Barker-Henderson (BH) perturbation theory and thermody¬ 
namic perturbation theory (TPT) for associating fluids in combination with BH perturbation theory are applied to 
the study of the structural properties and phase behaviour of the Derjaguin-Landau-Verwey-Overbeek (DLVO) 
model of lysozyme and y-cristalline aqueous electrolyte solutions. Predictions of the MSA for the structure 
factors are in good agreement with the corresponding computer simulation predictions. The agreement be¬ 
tween theoretical results for the liquid-gas phase diagram and the corresponding results of the experiment 
and computer simulation is less satisfactory, with predictions of the combined BH-TPT approach being the most 
accurate. 
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1. Introduction 

Among globular proteins, the equilibrium properties (structure factors, thermodynamic properties 
and phase behavior) of aqueus solutions of lysozyme and y-crystalline are the ones, perhaps, most thor¬ 
oughly studied. Numerous recent studies have been reported for these systems from both theoretical and 
experimental perspective (see (THg) and references therein). Much effort, in particular, has focused on 
the investigation of the phase behavior of the lysozyme solution. Experimentally, the liquid-gas phase 
coexistence for this system was described by Ishimoto and Tanaka 0. According to the later studies of 
Broide et al. (8), this phase coexistence is unstable with respect to crystallization. Worth mentioning is an 
important contribution due to George and Wilson O, who discovered the existence of the ‘crystallization 
slot’ for the values of the second osmotic virial coefficient in the vicinity of the liquid-gas critical point, 
where one might expect crystallization of the proteins. 

Since globular proteins can be viewed as colloidal macroions, most of the theoretical studies of pro¬ 
tein phase equilibrium are based on the concepts borrowed from the physics of colloids. As far as the 
interaction between two protein macromolecules is very complicated and, to a large degree, the known 
theoretical studies are based on the coarse-grained potential models (T0|. The simplest version of the 
model, the so-called one-component model, represents the effect of the solvent and electrolyte produced 
by the continuum approximation. Usually, the corresponding effective interaction is represented by the 
hard-sphere (or soft-sphere) interaction combined with long-range repulsive screened Coulomb interac¬ 
tion and short-range attractive van der Waals interaction. This is the model utilized in the Derjaguin- 
Landau-Verwey-Overbeek (DLVO) theory of the colloidal stability (11). Recent computer simulation stud¬ 
ies duniiaiii demonstrate that DLVO model, being not very accurate in reproducing the short-range 
behavior of the experimental structure factors at higher pH values, appears to be quite successful in 
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describing the phase behavior of the lysozyme solution. In particular, DLVO model was capable of repro¬ 
ducing the flat portion of the experimental phase diagram in the vicinity of the critical point. Further 
progress in the coarse-grained modelling of protein solutions is related to the introduction of the orienta- 
tionally dependent short-range attraction between colloidal particles 11411161 . The models of this type are 
aimed at a more detailed description of the protein molecules taking into account the existence of charged 
groups on their surfaces. Phase behavior and structural properties of these models have been studied in 
references Further steps in the detalization of the interaction in protein solutions are con¬ 

nected with the substitution of the one-component model with the one that contains multiple components 
and takes into account simple salt ions and solvent molecules. The presence of a simple electrolyte can 
be addressed within the framework of a highly asymmetric electrolyte solution model, which has been 
extensively studied using integral equation methods (see (23 and references therein) and by explicitly 
taking into account the effects of association 12611321 . Finally, the effects due to solvent molecules can be 
considered by extending the scheme developed for simple electrolytes and polyelectrolytes ISSMOl 
Despite a good performance of the DLVO model in predicting the phase behavior of the electrolyte 
solution of lysozyme (4l[3, the authors have not proposed any theoretical description for its equilibrium 
properties. The goal of the present study is to fill this gap and develop simple theories capable of describ¬ 
ing the structural properties and phase behavior of the model, at least at a qualitative level. With this 
goal In mind, we apply several simple hquid state theories, including the mean spherical approximation 
(MSA), second order version of the perturbation theory (TPT), for associating fluids in combination with 
BH perturbation theory, and critically assess their performance. The paper is organized as follows. In 
the section we formulate a potential model while in section we discuss the details of the MSA, BH 
perturbation theory and TPT. Our results are presented in section and our conclusions are given in 
section[5] 


2. The model 


DLVO model utilized in (4| [3 El) treats lysozyme solution as an effective one-component fluid of 
spherical particles with the number density p interacting via a pairwise additive potential which is only 
a function of the interparticle distance r. Pair interaction between particles consist of (i) a short-range 
attractive part of van der Waals term. 
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where (Tdlvo is an effective hard-core diameter and Ah is the Hamaker constant, and (ii) a Debye-Hiickel 
(DH) term 
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where Ze is a net charge of the lysozyme macromolecule in electronic units, eo is permittivity of vacuum, 
er = 86.765 - 0.3232 x TCC] E). The Debye-Hiickel screening length k Is defined by the expression 
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where I is the ionic strength of the solution, which takes into account the presence of ions due to the 
buffer and the added salts. In the present study, we neglect a weak dependence of DH potential on the 
temperature and assume for k and the values calculated at ambient conditions. 

The total DLVO potential is written as follows: 


i4)LV0(?') = 


oo. 


r < o'dlvo + ^. 

r S 0'DLVO + ^. 


(4) 


where the cut-off value S is introduced to avoid a singularity of the van der Waals contribution at r = 
ctdlvo end corresponds to the thickness of the Stern layer. 
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3. Theory 

The properties of the model at haitd are studied using MSA, BH perturbation theory and TPT for as¬ 
sociating fluids. MSA is a simple analytical approach, which is known for being capable of providing 
both structural and thermodynamic properties of a large number of model systems with sufficient ac¬ 
curacy. Relatively simple and possibly accurate description of thermodynamic properties of the model 
can be also achieved within the framework of the BH perturbation theory. Finally, to account for strong 
and short-range attraction between particles, which characterizes DLVO potential model we apply an 
appropriate combination of the BH perturbation theory and TPT for associating fluids. The accuracy of 
each of these approaches is evaluated via comparison of the theoretical predictions with the correspond¬ 
ing computer simulation and experimental predictions. 


3.1. MSA 

MSA consists of the Ornstein-Zernike (OZ) equation 

h[k) - c[k) + pc[k)h[k), (S) 

where c{k) and h[k) are Fourier transforms of the direct and total correlation functions c(r) and h(r), 
respectively, and MSA closure relation 


f cir] =-p[V^awir) + VuH(r)], r^aoLvo + S, , . 

( h(r) = -l, r<( 7 DLvo + ^- 

For the sake of analytical description, we approximate the vdW part of the potential Vvdwf'') using 
one-Yukawa potential, i.e., 


1/ r, w,, A exp[-KY(r-((TDLVO + '5)] 

KrdwC'') = Wfr) = - Ay-. (7) 

r 

As a result, the pair potential outside the hard core is represented by the two-Yukawa potential and 
we have used an analytical solution of the MSA due to Blum and Hoye 141]. The structural and ther¬ 
modynamic properties have been calculated utilizing this solution and using a closed form of analytical 
expressions presented in 1^ . 

3.2. Barker-Henderson perturbation theory 

Here, we utilize the second-order BH perturbation theory. Within the framework of BH perturbation 
theory 03, Helmholtz free energy F of the system per particle is given by the following expression 

p jcjideal pHS p_^ p^ 

NkjiT ~ Nk^T Nk^T Nk^T Nk^T ’ 
where fcs is the Boltzmann constant, is the ideal gas free energy 

j^ideal 

-= ln(Ap]-l, (9) 

NkjiT w - 

A is the thermal de Broglie wavelength, is the hard-sphere Helmholtz free energy and Fi and F 2 are 
the first- and second-order contributions to Helmholtz free energy of the system. Here, for we have 
used the Carnahan-Starling extension 

F^s _ 3 ,j2 

-= — -V - (10) 

NkTiT (1-77)2 
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where ?? = f pc^Q, (Jq = ctolvo + ^ and for Fi and F 2 , we have 

00 
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( 11 ) 


( 12 ) 


where is the compressibility of the hard-sphere fluid. Here, we use the corresponding Carnahan- 
Starhng expression l44l . 

Using the above expression § for Helmholtz free energy of the system, all the rest thermodynamical 
properties (pressure and chemical potential) can be derived using the standard thermodynamical rela¬ 
tions. As for any perturbation theory, the BH model works best if the attractive potential is not large in 
magnitude. For potentials with strong attractions, other approaches are more suitable. 


3.3. Thermodynamic perturbation theory for associating fluid 

To account for the strong attraction seen in the DLVO potential of lysozyme, we combine the BH 
perturbation theory and thermodynamic perturbation theory (TPT) for associating fluids with spherically 
symmetric interaction (45). Following the earlier studies (26l|29]|46), the total pair potential of the system 
is represented as a sum of the reference and associating pieces, i.e., 

VoLVO(r) = <LVO+ ^DLVO> (13) 
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Here, FbLvo(o'o) < ^ assume that is temperature dependent. A particular choice for 

the potential splitting parameter is discussed below. According to the splitting of the total DLVO 
potential l |T^ , the system Helmholtz free energy within the framework of the TPT (43 is as follows: 

(IS) 


where for the free energy of the reference system F’’®^ we have used the second-order expression ^ with 
DLVO pair potential Volvo (^) substituted by the reference potential V^'j^plr), and thus for the associative 
part F*®®' we have 
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m is a maximum number of bonds per particle allowed and is the contact value of the radial 

distribution function of the reference system. The latter quantity was obtained using the parametrization 
of the contact value of the radial distribution function for the hard-sphere square-well fluid (47). 

The knowledge of the free energy of a system enables one to calculate thermodynamical prop¬ 
erties of interest using standard thermodynamical relations. The version of the TPT approach discussed 
above has two input parameters, i.e., the maximum number of bonds m and the potential splitting pa¬ 
rameter . In general, the number of nearest neighbours for the model at hand could be up to 12, thus 
one can assume m - 12. In this case, the value of which ensures the saturation of the associative 
potential Vqlvo(c) for m - 12, can be chosen to be relatively large in comparison with the contact value 
of the DLVO potential Vdlvo(o’o)- However, according to the earlier studies [4^, ‘single bond’ approxima¬ 
tion utilized here is accurate only for relatively small values of m. In addition, the probability of bonding 
of 12 particles simultaneously at the densities where the liquid-gas separation occurs is small; thus, the 
optimal choice for m and requires a certain compromise. In the present study we assume that m-3 
and the reduced value of the associative potential at the contact T is constant. 

In addition, to provide an accurate description of the reference system at the critical point, we also as¬ 
sume that Vq®®' = fcs Ter, where Ter is the critical temperature. Combining the latter two assumptions, we 
have 


T r(as) 
^ _ 

fcsT 


= HdlVO(O'o) 


1 


Tct 


- 1 . 


( 20 ) 


For T - Ter we have = - fcs ?cr. i-e., for the critical temperature, the minimum value of the reference 
system potential Is equal to -k^ T^r- For this value of the potential minimum, the BH approach 

is expected to be sufficiently accurate. Although our choice for m and Vq®®’ is rather empirical, the the¬ 
ory proposed Is self-contained because there is no need in the input from outside. In particular, critical 
temperature T^r (and critical density pcr)> which enter the expression for Vg®®' 1 20 1 , is obtained as usual 
from the solution of the set of two equations, which requires the first and the second derivatives of the 
pressure with respect to the density to be equal to zero, I.e., 


dP[p, Ter) I 

dp lp=pcr 
d^Pip, Ter) I 

dp2 lp=pcr 


( 21 ) 


4. Results and discussion 

In this section we present our numerical results for the structural properties and phase behavior of 
lysozyme and y-crystalline aqueous electrolyte solutions. In both cases, we use the DLVO model with 
o’DLVO = 37.08 A for lysozyme and ctdlvo = 37.8 A for y-crystaUine. For the Stern-layer thickness and 
Hamaker constant, the following values have been used O: 5 = 1.8 A and Ah = 18.8 kj/mol. 

The structural properties of the lysozyme solution [static structure factor S(k] and radial distribu¬ 
tion function g(r)] were calculated using MSA |^. Parameters for the Yukawa potential |7} representing 
Vvdw(?‘) were chosen using ‘best eye fit’ method supplemented by the equality of both potentials at the 
contact distance, i.e., Vy{(Jo) - 14 dw(o’o)- In figure[^we compare VyW and Vvdw(^) with the following 
choice of the Yukawa potential parameters: Ay = 375.19 kJ/mol-A and ky = 0.45 A“^. In figure|^we com¬ 
pare the corresponding theoretical and computer simulation O results for the structure factor S(fc) at 
the values of pH used in the corresponding experimental study, i.e., pH = 2.8, 4.2, S.07. These values 
of pH correspond to the following three values of the ionic strength (see table 1 of reference O), I.e., 
I - 0.035 M, 0.081 M, 0.102 M, respectively. In addition, we also show the corresponding experimental 
results m for S(fc). Very good agreement is observed between theoretical and computer simulation re¬ 
sults for the structure factor at all values of pH and I studied. Although DLVO model is not very accurate 
in describing the short-range behavior of the structure factor at higher values of pH (figurej^, its predic¬ 
tions for the phase behavior are in reasonable agreement with experimental predictions (see figure [^. 
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Figure 1. Van der Waals contribution to the DLVO potential (line with triangles) [see j^] and Yukawa 
substitution (full line) of van der Waals term. 


These features of the DLVO model 0 are known and have been discussed earlier in figure we 
show our theoretical results for the liquid-gas phase diagram of the lysozyme solution at pH = 6.0 and 
ionic strength / = 0.6 M (black color) and y-crystalline aqueous electrolyte solution at pH = 7.1 (which 




Figure 2. Structure factor of the lysozyme aqueous electrolyte solution at different pH. 
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Figure 3. (Color online) Liquid-gas phase diagrams of the lysozyme aqueous electrolyte solution at 
pH = 6.0 and ionic strength 7 = 0.6 M (black color) and y-crystalline aqueous electrolyte solution at 
pH = 7.1 and I = 0.24 M (red color). Experimental values for the critical concentration and tempera¬ 
ture for lysozyme (49j and y-crystaUine solutions (S^ are (230 + 10 mg/ml, 273 K) and (289 + 20 mg/ml, 
311 K), respectively. 


correspond to the protein overall charge Z - 1) and ionic strength I = 0.24 M (red color). These results 
were obtained using second-order BH theory and our version of the TPT (T^ . Unfortunately, for the 
employed parameters of the DLVO model, the MSA does not have a convergent solution in the range of 
temperatures and densities, where one would expect the location of binodals. In addition, in the same 
figure we show computer simulation results and the results of the experiment (4^[^. Here, computer 
simulation predictions are much more accurate in comparison with the predictions for the structure 
factors demonstrated in the previous figure. Our results obtained using the BH perturbation theory, are 
much less accurate, giving the values for both critical temperature and critical density that are too large. 
In addition, BH phase diagram is too narrow in comparison with the computer simulation and experi¬ 
mental phase diagrams. Predictions of our combined BH and thermodynamic perturbation theory are in 
reasonable agreement with the computer simulation data. Good agreement is observed for the critical 
temperature. Slightly less accurate are the predictions for the critical concentration. The overall shape of 
the phase diagram is still too narrow, although it is now much closer to the shape of computer simulation 
phase diagram. 


5. Conclusions 

In this article we studied the structural properties and phase behaviour of the DLVO model of lyso¬ 
zyme and y-crystalline aqueous electrolyte solutions using MSA, second-order BH perturbation theory 
and a combined approach based on the BH theory and TPT for associating fluids. Theoretical results are 
compared with computer simulations and experimental results. Predictions of the MSA for the structure 
factor of lysozyme solution are in good agreement with the corresponding computer simulation predic¬ 
tions. However, MSA does not have a convergent solution in the range of the temperatures and densities, 
where one would expect the location of the corresponding MSA binodals. We conclude that MSA is in¬ 
appropriate for the phase behavior of the DLVO-type models of lyzosyme and, perhaps, other globular 
proteins. Among the theories used to describe the phase behaviour of the lysozyme and y-crystalline 
solutions, only a combined BH-TPT approach provides a reasonable qualitative agreement with com¬ 
puter simulations and experimental description. We expect that a further improvement of the theory 
can be achieved using association concepts in combination with a more detailed description of the pro- 
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tein molecules taking into account the existence of the charged groups on their surfaces. Corresponding 
studies are underway and results will be reported in due course 1^ . 
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TeopeTi/iMHe AOC-niA>KeHHji (|)a 30 BoT noeeAiHKi/i MOAe/ii 
BOAHi/ix posMi/iHiB e/ieKTpo/iiTlB /lisoi^i/iMy i raMMa-Kpi/iCTa/imy 

P. MehbHi/iK 

iHCTntyT c()i3i/iKn KonflencoBaHi/ix ci/icreM HaqiOHaabHoT aKafleMli HayK VKpaTHM, 

Byn. I. CBeHqiqbKoro, 1,79011 /IbBlB, VKpaiHa 

CepeflHbOCcJepi/iHHe naOai/ixeHHB (CCH), reopin sOypeHb flpyroro nopn^Ky BapKepa-feHflepcoHa (Bf) I repMO- 
AMHawlHHa reopia sOypenb (TT3) flan acoqlaiMBHi/ix pifli/iH b KOMOiHaqii 3 leopieK) sBypenb Bf 3acTOCOBaHi 
p,o BMBHeHHB CTpyKTypHMx BaacTi/iBOCTePI ra (jjaaoBOi noeefllHKn MOfleai /lepnriHa-/laHAay-BepBeB-OBep6eKa 
(fl/lBO) BOflHMX po3HMHiB eaeKipoahlB aiaopnMy I raMMa-Kpncraamy. PeayabTarn CCH flan crpyKiypHnx (fia- 
KToplB flo6pe yaroflxyKDTbca 3 BlAnoBlflHiaMn KOMn'rorepHiaMn poapaxyHKaMM. y3rofl>KeHHa Mi>K TeoperianHiaMn 
peayabTaraMM flan (jjaaoBOi fliarpaMia piflMHa-raa I BiflnoBiflHiaMn eKcnepnMeHiaMn ra KOMn'rorepHiaMn cnr/iy- 
aaqiaMM e Meniu aaflOBlabHi/iM. Kpautoro yaroflxeHHa flan (fiaaoBnx fliarpaM fl03B0aae flOcarHyri/i KOwSiHOBa- 
huPI Bf-TTS niflxifl. 

iOiioMOBi cnoBa: MOfle/ib/I/IBO, noreHLiia/i HDKasM, BOflnipoa^nHi^/liaoi^i/iMy i raMMa-Kpi/icran'iHy, 
cepeflHbOdpepi/mHe Ha6ni/i>KeHHa 
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